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ABSTRACT 

Context. Unlike gas-phase reactions, chemical reactions taking place on interstellar dust grain 
surfaces cannot always be modeled by rate equations. Due to the small grain sizes and low flux, 
these reactions may exhibit large fluctuations and thus require stochastic methods such as the 
moment equations. 

Aims. We evaluate the formation rates of H 2 , HD and D 2 molecules on dust grain surfaces and 
their abundances in the gas phase under interstellar conditions. 

Methods. We incorporate the moment equations into the Meudon PDR code and compare the 
results with those obtained from the rate equations. 

Results. We find that within the experimental constraints on the energy barriers for diffusion and 
desorption and for the density of adsorption sites on the grain surface, H 2 , HD and D 2 molecules 
can be formed efficiently on dust grains. 

Conclusions. Under a broad range of conditions, the moment equation results coincide with 
those obtained from the rate equations. However, in a range of relatively high grain temperatures, 
there are significant deviations. In this range, the rate equations fail while the moment equations 
provide accurate results. The incorporation of the moment equations into the PDR code can be 
extended to other reactions taking place on grain surfaces. 

Key words. Gas-grain models - deuteration processes - molecular hydrogen formation - mo- 
ment equations 



1. Introduction 

Most of the detected interstellar molecules are in the gas phase. The chemical balance of gas phase 
species in interstellar clouds is commonly described by a set of rate equations, using the known 
rate coefficients of gas-phase reactions. Many studies have been devoted to the resolution of the 
chemical rate equations to derive the compositions of interstellar clouds within specific physical 
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conditions. Two main approaches are used by modellers: the time-dependent approach and the 
steady-state approach. In the time-dependent approach, one follows the time evolution of the abun- 
dances of various atoms and molecules by numerically integrating the corresponding set of first 
order coupled differential equations, from specific initial conditions. In the steady state approach, 
one directly solves the set of algebraic equations which is obtained when the time derivatives of 
the abundances are set to zero. Comparison with observations allows then to test the relevance of 
the various hypotheses. However, some important reactions which give rise to the formation of 
molecular hydrogen, ice mantles and certain organic molecules do not take place in the gas phase, 
but on the surfaces of dust grains. To account for these surface processes, one needs to incorporate 
the grain-surface reactions into the gas phase models of interstellar chemistry. It is thus attractive to 
use rate equations for the surface reactions, which enable to couple surface reactions to gas-phase 
reactions in a straightforward way. This approach is indeed valid in the case of sufficiently large 
grains, when the number of reactive atoms and molecules of each species on a grain is sizeable. 
However, in the limit of small grains, when the number of reactive atoms and molecules on a grain 
is small, the rate equations are not always suitable. Indeed, rate equations simply ignore the fluctua- 
tions in the number of reactive species on the grains. This problem was discussed by several authors 
over many years ( Tielens & Hagen|1982 Charnley et al.|1997 Casein" et al.|19 98 ; Shalabiea et al. 



1998 Stantcheva et al. 2001 2002} . To overcome these difficulties, a master equation approach has 



been proposed and was specifically applied to the formation of molecular hydrogen. This approach 
is suitable for the simulation of diffusive chemical reactions on interstellar grains in the accretion 



limit ( Biham et al. 2001 Green et al. 2001 Biham & Lipshtat 2002 ). The master equation approach 
takes into account the discrete nature of reactive species on the surfaces of grains as well as the fluc- 
tuations in their populations. This approach was later used for a larger chemical network on grains, 



leading to the formation of methanol and its deuterated versions (Stantcheva et al. 2002 Stantcheva 



& Herbst 2003 ). However, the number of equations included in the master equation increases expo- 



nentially with the number of species that are reactive on grain surfaces. Therefore, it is not suitable 
for incorporation in codes that include complex networks of grain-surface reactions. In addition to 



the master equation method, a Monte Carlo method was proposed by Tielens & Hagen ( 1982 1 and 



applied to interstellar modelling by Charnley (2001 ), Caselli et al. (2002 1 and Stantcheva & Herbst 



(2003). Recently, Cuppen & Herbst (2005) developed random walk models for molecular hydro- 
gen formation on grains that take into account the effect of surface roughness on the diffusion and 
reaction rates. They simulated these models using Monte Carlo methods and showed that surface 
roughness tends to broaden the range of temperatures in which molecular hydrogen formation is 
efficient. The advantage of the random walk models is that they enable to account in more detail 
for the microphysics of the grain surface. In particular, they enable to take into account the whole 
distribution of binding energies and diffusion barriers for H atoms. Detailed models of this type 
provide useful insight but it is not feasible to include them in large models of interstellar chemistry 
due to their complexity. In contrast, the rate equation and moment equation methods used in this 
paper, include a single energy barrier for each process. Finally, a semi-empirical approach, known 



as the modified rates method, was proposed by Caselli et al. ( 1998). This method is easy to em- 



ploy and has been used with mixed success by |Shalabiea et al.| ( |1998[ ), [Stantcheva et al.| ( |2002[ ) and 
Caselli et al. (2002). Garrod (2008 ) has studied a different version of the modified rate method. 
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The recently proposed moment equations method provides efficient stochastic simulations of 
grain-surface chemistry (Barzel & Biham 2007a b). The method consists of only one equation 
for the average population size of each reactive species on a grain and one equation for the rate 
of each reaction. It consists of a set of coupled ordinary differential equations, that resemble the 
rate equations. Therefore, they can be easily coupled to the rate equations of gas-phase chemistry. 
Moreover, unlike the rate equations, the moment equations are linear, and thus easier to handle 
both in the time-dependent and in the steady state approaches. The moment equations were tested 
for the reaction network that gives rise to ice mantles on grains, that consist of water ice, carbon 
dioxide and methanol (Barzel & Biham 2007a b ). The stability properties of the moment equations 
as well as the accuracy of their steady-state solution under astrophysically relevant conditions were 
examined by extensive computer simulations and by comparison to the master equation results. 



In this paper we incorporate the moment equations into the Meudon PDR code (Le Petit et al. 



2002 2006[ l, a stationary model of Photon Dominated Regions. The model considers a stationary 



plane-parallel slab of gas and dust illuminated by an ultraviolet radiation field coming from one 
or both sides of the cloud. It solves, in an iterative way, at each point in the cloud, the radiative 
transfer in the UV, taking into account the absorption in the continuum by dust and in discrete 
transitions of H and H2. Explicit treatment is performed for C and S photoionization, H2 and HD 
photodissociation, as well as CO (and its isotopomeres) predissociating lines. The model also com- 
putes the thermal balance, taking into account heating processes such as the photoelectric effect 
on dust, cosmic rays and chemistry. It also accounts for the cooling due to infrared and millimeter 
emission of the abundant ions, atoms, and molecules. The chemistry is solved under steady state 
conditions and the abundance of each species is computed at each point. The excitation states of 
a few important species are then computed. The column densities of these chemical species and 
their emissivities/intensities are then calculated. To examine the applicability and relevance of the 
moment equations within the PDR code, we consider a simple network of grain surface chemistry, 
involving only H, D, H2, HD and D2. We compare the results obtained from the PDR code with 
the moment equations with those obtained when the rate equations are incorporated in the same 
code. Unlike previous studies of molecular hydrogen formation that assumed a single grain size, 
the moment equations enable us to take into account the full distribution of grain sizes. We show 
that in case of relatively high grain temperatures, the rate equations are not suitable and the moment 
equations should be used. We also demonstrate the importance of the Langmuir rejection effect at 
very low grain temperatures, where the grain surfaces are saturated by hydrogen atoms and reaction 
rates are low. 

The paper is organized as follows. In section[2]we briefly review the gas-phase and surface reac- 
tions involved in the formation of H2, HD and D2 molecules. We present in section[3]the moment 
equations for the H-D system and describe their incorporation into the PDR model. Simulation 
results are displayed in section|4] The results are summarized and discussed in section[5] 



2. Formation of Molecular Hydrogen and its Deuterated Versions 

2.1. Astrophysical context 

No efficient gas-phase mechanism is at work to form molecular hydrogen in the gas phase. 
Therefore, reactions on grain surfaces are needed in order to account for the observed abundance 
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of H2.[Hollenbach & S alpeter] ( | 1 97 1 ) l have been among the first to quantitatively describe the for- 
mation of molecular hydrogen on the surfaces of spherical dust grains. Here we assume that the 
grain size distribution follows a power law of the form 

n(a) = ca- q , (1) 

as suggested by Mat his, Rumpl, & Nordsieck| ( |1977| l. The prefactor c is given by 
3 1.4-m H -G-(4-g) 

c ; ~, «h (2) 

4;r p(a 4 -i _ fl 4 -<?) 

AA"max "min^ 

where p (g cm~ 3 ) stands for the volumic mass of the grains, while a max and a m ; n are the upper and 
lower cutoffs of the grain radii distribution, respectively. We denote the proton density in the gas 
phase by «h (cut 3 ) and the dust to gas mass ratio by G = 0.01. In this case, one can formulate the 
formation rate of H2 as a function of various grain parameters such as the dust to gas mass ratio, 
the minimum and maximum values of the grain radii and the volumic mass of the grains. If the 



exponent of the power law is q = 3.5, analytical formulae are obtained as derived by Le Bourlot 



et al. ( 1995] >. Le Petit et al. (2002) have specifically studied the D/HD transition occuring at the edge 



of a translucent cloud and in dense photon dominated regions. They have also studied the formation 



of HD molecules on the surfaces of dust grains, following the procedure of Le Bourlot et al. ( 1995 ). 
The resulting formulae are very similar to the pure molecular hydrogen case. Following previous 
considerations by |Watson| ([1974} and |Dalgarno, Black, & Weisheit| ( |1973| l, |Le Petit et al] ( |2002| l 



found that HD is much more efficiently formed in the gas phase, in a succession of reactions 
initiated by cosmic ray ionization of atomic hydrogen. When protons are formed, the sequence of 
reactions involved is the following: 

H + + D ^ D + + H. (3) 

The forward reaction is exothermic by an energy of 3.7 meV, corresponding to 43 K, so that the 
reverse reaction may take place at moderate temperatures, reducing the D + formation efficiency. 
D + reacts efficiently with H2 to form HD: 

D + +H 2 -> HD + H + . (4) 

The reverse reaction may also occur, but the endothermic barrier is now 40 meV, equivalent to 
464 K, so that the corresponding probability is negligible at low temperatures. The main destruction 
channel of HD in diffuse environments is photodissociation, which has been computed by |Le Petlt| 
et al. ( 2002| l. This simple gas-phase scheme accounts for the observed HD/H2 abundance ratio in 



diffuse clouds and translucent regions, as found from Copernicus and FUSE observations (Watson 
1974 Dalgarno et al. 1973 Lacour et al. 2005 ). In dark and cold regions, HD is the major repository 



of deuterium. It reacts with H3 to form H2D + which is the starting point of a rich chemistry of 
deuterium fractionation ( |Millar et al.p989) |Roueff et al.pOOOl |Caselli et al.]2002l |Stantcheva &| 
|Herbst|[2003| |Roberts et al^[2003l |Roueff et al.pOOB) . Although the formation of HD in diffuse 
clouds occurs primarily in the gas, it was recently proposed that under certain conditions, grain- 
surface reactions may also contribute to an enhanced production of HD and D2 molecules. The 
proposed mechanism is based on the assumption that D atoms stick more strongly than H atoms so 
that their desorption rate is lower. Such an isotope effect has been observed in various experimental 
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5 



situations (Koe hler et al.|1988||Ho ogers 1995; Amiau det al.|2007) . As a result, the residence time 
of D atoms on grains is expected to be longer than that of H atoms. Then, the D/H abundance ratio 
on the grains is enhanced compared with the gas-phase ratio and newly adsorbed H (or D) atoms 
are more likely to find D atoms already residing on the grains. This may give rise to an enhanced 
production of HD molecules. 



2.2. The interaction between hydrogen/deuterium atoms and dust grains 

In order to quantitatively describe the formation of H2, HD and D2 molecules on grain surfaces, one 
has to introduce several hypotheses and parameters, which we recall below. The typical velocities 
of H and D atoms in the gas phase are given by vu and vq, respectively. Under the simplifying 
assumption that all grains are of the same radius, a, one can compute the numerical density of 
spherical dust grains n gr (cm -3 ) as a function of riu- 

3 • m H • G 

'V = 1.4- r « H - (5) 

4-n ■ a' ■ p 

For a power law distribution of grain sizes with the exponent q — 3.5 (Ma this et al.||1977 l, the 
numerical density of the grains becomes: 

3 • ntii ■ G ■ (a" - a 2 5 ) 

1 a n \ max mm 7 /r\ 

n„ r — 1.4 • =-z • «H, (o) 

207T ■ p ■ (aminflinax) ' ( V fl max ~ V a min) 

The number of adsorption sites on a grain is denoted by S . Their density on the grain surface, s 
(sites cm -2 ), is given by s — S /4na 2 . The density of adsorption sites, s, and the distance between 
adjacent sites, d, are related through d 2 = l/s. The fluxes F H and F D (atoms s" 1 ) of H and D atoms 
onto the surface of a single grain are given by Fi = n\V\cr, where I = H, D, respectively and cr is 
the cross section of a grain, namely cr = na 2 . 

The atoms stick to the surface and hop as random walkers until they either desorb or recombine 
into molecules. The desorption rates of H and D atoms on the surface are given by 



Wi = vi • exp 



kTgr j 



(7) 



where v\ is the attempt rate, £^ es is the energy needed for desorption of an atom of isotope I and 
T gr is the grain temperature. The hopping rates of the atoms on the surface are 



«i = vi • exp 



(8) 



where Ej lff is the activation energy barrier for diffusion of the isotope I. The rate a\ is the inverse 
of the residence time ti of an atom of isotope I in a single adsorption site. The sweeping rate 
A\ = ai/S is approximately the inverse of the time Sti required for an atom of isotope I to visit 
nearly all the adsorption sites on the grain surface. Since the D atom is twice as heavy as the H 
atom, its ground state energy within an adsorption site on the surface is lower. Due to this isotope 
effect, the desorption barrier for D atoms is assumed to be higher by 5 meV than for H atoms. For 
the diffusion barriers, we assume £'ji ) lff = Ej? ff because the diffusion barrier balances the zero point 
energy of the potential well with that of the saddle point (or transition state), while desorption does 
not possess such saddle point. Note that in the analysis of the experimental results, presented in 
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Table 1. Parameters for the interaction of hydrogen and deuterium atoms with dust grains of dif- 
ferent compositions and surface morphologies. Numbers in parentheses refer to powers of ten. 



Parameter 


Amorphous carbon 
Katz (1999) 


Olivine 
Katz (1999) 


Amorphous silicate 
Perets (2007) 


Low density ice 
Perets (2005) 


P (g crrr 3 ) 
s (crrr 2 ) 
d(A) 

E^ s (meV) 
E£ s (meV) 
Ef s (meV) 
v\ 


2.16 

5(13) 

14.14 

56.7 

61.7 

44.0 

10 12 


3" 

2(14) 

7.07 

32.1 

37.1 

24.7 

10 12 


3.5 

7(14) 

3.78 

44.0 

49.0 

35.0 

10 12 


0.94 

5(13)" 

14.14" 

52.3 

57.3 

44.5 

10 12 



a Assumed sensible values, as they are not given in the cited papers. 



|Katz et aT] ( | 1999[ > and |Perets et al.| ( |2005| ), no distinction was made between H and D atoms, in 
order to keep the number of fitting parameters small. 

The possibility of tunneling of the adsorbed atoms between adsorption sites has recently been 
studied by |Cazaux & Tielens| ( [2004l ). 

The tunneling probability depends on the distance between adjacent sites and on the energy 
barrier, which is taken as the diffusion barrier. For the sake of completeness, we present the corre- 
sponding tunneling rate, which is given by (Eisbgerg 1961 1 : 



1 + 



sinh 2 (d/A DB ) 



4E (Ef iff - E) 

where A^b corresponds to the De Broglie wavelength and is given by 



(9) 



tDB 



(10) 



E) 



and E is the kinetic energy of the adsorbed atom, determined by the grain temperature, namely 
E = k B T gr . The efficiency of the tunnelling process is sizeable only when the ratio d/A DB is of order 
1 or less. There is no specific reason to assume that the adsorption sites are uniformly distributed. 
However, for simplicity we assume that the distance d between adjacent sites is fixed and given 
by d 2 = l/s, as explained above. We have evaluated the tunneling rate for the parameters used in 
this paper, and it was found to be negligible. Therefore, in the simulation results presented in this 
paper, the mobility of H and D atoms on grains is only due to thermal hopping, while tunneling is 
ignored. 

In Table[T]we display the parameters that describe the interaction of H and D atoms with grains 
of different compositions and surface morphologies, as reported in the literature. These parameters 
are based on a series of experiments and subsequent analysis, reported in Katz et alT]( |1999) l, |Perets| 
and '. 



et al. 



(2005 



Perets et al. 



12 ,-1 



(2007 1. The attempt rate vi is often assumed to be equal to 10 s 



It may also be derived from an harmonic oscillator model, where vi = J2E^Jm j{n d). 
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The reported values of the desorption and diffusion barriers put severe constraints on the grain 
temperature range over which molecular hydrogen formation may occur. If one considers a distri- 
bution of grains of various sizes, such as the power law derived by Mathis e t al.| ( |1977[ ), one has 
to perform the average of the various quantities over the size distribution, keeping in mind that 
different temperatures may pertain to grains of different sizes. In the present paper, we assume that 
the grain temperature depends on the surrounding radiation field but not on the grain size. Under 
this assumption, in a given location in the cloud, grains of all sizes have the same temperature. 
In the simulations reported in this paper we use the parameters of amorphous carbon. Amorphous 
carbon is assumed to be a primary component of interstellar dust. Its surface properties are thus 
suitable for the simulation of molecular hydrogen formation in relatively warm regions in which 
the grain surfaces are not covered by ice mantles. It is also observed that the parameter values 



of amorphous carbon are close to those of low density amorphous ice, reported by Perets et aL] 
( p005] >. Therefore, we use the same set of parameters for the interaction of hydrogen atoms with 
grain surfaces, independently of whether these grains are expected to be bare or covered by ice 
mantles. 



3. Rate equations, master equation and moment equations 

3.1. Rate equations 

After being adsorbed onto the surface, the H and D atoms hop between adsorption sites until they 
either desorb or form new molecules. The number of atoms of isotope I on the surface of a grain is 
denoted by iVi : . The rate equations account for the expectation values (Ni), of the population size 
of isotope I on a grain of a given radius where I = H or D. 



3.1 .1 . Basic processes 

In the present version of the PDR code, we introduce three main physical processes to describe 
the formation of H2, HD and D2 on grain surfaces: adsorption, diffusion-mediated reaction and 
desorption: 

- Adsorption 

H + grain -» H: (fc^ins -1 ) 
D + grain -> D: (k^ in s _1 ), 

where H: and D: stand for the hydrogen and deuterium atoms which are adsorbed on 
the grain surface. The corresponding rate equations are: 

^ = +*5.-CH) (id 
^ = (12) 

The adsorption coefficients k^ ds and k^ are directly proportional to the grain cross section 
<x = na 2 and to the sticking probability y, namely 

*5i = ^)"g r = 7 v ^a 2 n gr . (13) 
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For a distribution of grain sizes, one has to perform the average on crn gr over the distribution 
function, which leads to the following expression: 



,h , F u x , „ 3 • m H • G 

n(H) 4n ■ p ^a min a„ 



Cs=<I7^"gr>=y^H-1.4- _ " _ -n H (14) 



Similar equations hold for deuterium adsorption. 
Reaction on the grain surface due to diffusion 

The diffusion-mediated formation of molecular hydrogen and its deuterated versions is de- 
scribed by: 

(k"; f in cnrV 1 ) 
(k™ in cnrV 1 ) 
(k^incm-V 1 ). 

Assuming that all formed molecules are directly released into the gas phase, the formation rates 
of H2, HD and D2 are given by : 

= kf^n(H :) 2 = A H (N H .)\ r (15) 
= Cf»( H On(D :) = (A H + A D )(N H )(N D .)n gr (16) 



H : 


+ 


H : 




H 2 


H : 


+ 


D : 




HD 


D : 


+ 


D : 


— > 


D 2 



d -^=k^n(D :f = A D (N D .)\, 



= k*n(D:) 2 =A D (N D .)\ r , (17) 



where (Nh) ({Nd-}) is the number of adsorbed hydrogen (deuterium) atoms on a single grain. 
The rate coefficients corresponding to these surface reactions can be obtained from the rela- 
tions: n(H :) = (A^H:> ■ n gr and n(D :) = (No) • n gr . These rate coefficients take the form 



, h 2 _ A H (N H y _ Ah 
surf " n(H :)2 Hgr ~ n gr 

(A H +A D )(N H: )(N D: ) A H +A E 



k ™ " n(H :)n(D :) " gr ~ n gr 
, d 2 _ A D <jV D: ) 2 _ Ad 
surf " «(D:)2 Hgr ~ n gr - 

These formulae are valid only if the grains have the same temperature because the diffusion 
coefficient is temperature dependent. 

Desorption 

The desorption processes 
H: -> H + grain (k^ins -1 ) 

D: -» D + grain (kJjMns- 1 ), 
are described by the equations 



dn(U:) =-C s n(H0 (18) 



dt 



(19) 



where 
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**, = w h (20) 
**, = Wb- (21) 

The overall variations of adsorbed H and D atoms per cm 3 in the rate equation formalism are: 



= k» s n(H) - W»n(H :) - 2^„(H :f - ^±^„(H :)n(D :) (22) 
dn(T> :) D A D 2 (A H +A D ) 

— 3: — = k *dA D ) - ^d«(D :) - 2 — n(D :) n(H :)n(D :) (23) 

(XX Y\gy TXgf 

3.1.2. Rejection effects 

As the number of adsorption sites on a grain is finite, one has to take into account the possibility 
that hydrogen atoms that impinge upon already occupied sites may be rejected rather than adsorbed 
on the grain. This mechanism is often referred to as Langmuir rejection. When the rejection is taken 
into account, the flux terms are modified according to 

Ff = F,(l - (24) 

where I stands for H or D. Since the abundance of D atoms is very small (1.5 1CT 5 ) compared to 
H atoms, it is sensible to consider first, for simplicity, the case in which the rejection is only due to 
already adsorbed H atoms. 

- Rejection only due to adsorbed H atoms 

We first consider the case of rejection only due to adsorbed hydrogen atoms. The rate equation 
for adsorbed H on a single grain becomes: 

= F H - (W H + y 1 ) (N H ) - 2A H <jV H :> 2 - (A H + A D )(^ a ><^D:>. (25) 

The net result is an apparent increase of the desorption term by /h = F^/S . When the incoming 
flux of deuterium atoms is similarly modified, the rate equation for the adsorbed deuterium 
becomes: 

" <M:) ~ F D - (w H + ^^$4 W> - 2A D <jV D: > 2 - {A H +A D )(N H .)(N D ..). (26) 



dt \ S (N D ;) 

The correction due to rejection introduces a term proportional to (Nh:)/{Nd-). The correspond 
ing equations for the quantities expressed in volumic density become the following: 

dn(H:) H ( yv H n(H)\ A H 2 (A H +A D ) 

= K<te n W -\ w n + — l«(H :) - 2— n(H :) «(H :)«(D :) (27) 



dt aas \ 4s 



igf fig r 



dn(D:) D I yu D n(D) n(H :) \ A D 2 (A H + A D ) 

= ^ads"( D )-| w d + — — |"(D :)-2— n(D :) : n(H :)n(D :)(28) 



dt aas \ 4s n(D :) / n gr n gr 

- Rejection by adsorbed H and D atoms 



We now include the adsorbed deuterium atoms in the rejection term according to Eq. (24i. As 
the equations are very similar, we only display the resulting rate equations: 

dn(U:) H ff A H 2 (^h+A d ) 

— 3: — = ^ads"( H ) _ w h "( H : ) - 2 — n(H:) n(H :)n(D :) (29) 



10 
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dn(D :) 
dt 



k^n(D) - Wfn(D :) - 2— n(D :) 



D 



(A H + A D ) 



"gr n gr 

The values of the effective desorption coefficients are : 



«(H :)n(D :). 



W eff = Wh + Wh"(H) + yv H n(H) n(D 



4s 



4s n(H 



= W D + 



yu D n(D) yvDti(D) n(H 



4s 4s n(D : 

The production rates of H2, HD and D2 (in units of cm _3 s _1 ) are: 



(30) 

(31) 
(32) 



R Hl = (AhAvMH :) 2 

Rhd = ((A H +A D )/n ?r )n(H:)n(D :) 

Rd 2 = (A D /n gr )n(D :f . 



This is a generalized form of the equations given in Lipshtat et al. (2004). 



(33) 
(34) 
(35) 



3.2. Master and moment equations when only the rejection term for H atoms is included 
3.2.1. Master equation 

When the number of adsorbed atoms is small, the rate equations may become inaccurate. In this 
case one may consider using the master equation or the moment equations derived from it. The 
master equation describes the temporal evolution of the probabilities P(Nh: , No.) that Nh : hydrogen 
atoms and Nd : deuterium atoms reside on the surface of a given grain. Here we write the master 
equation when only the Langmuir rejection due to hydrogen atoms is taken into account. It takes 
the form 



dP(N H ;,N D ;) 

dt 



= F H 
+ F D 



JP(A?H: - l,JV D: )-(l - ^p)p(jV H: ,jV D: ) 



(l - ^PiNu.Nu: - l)-(l - ^)P(^V H :,^V D: 

+ W H [(N H: + l)P(N H . + 1,N D ..) - N h P(Nh:,N D: )] 

+ W D [(iVb: + 1)P(^Vh:,^V D: + 1) - JV D:J P(JV H :, AW] 

+ A H [(Nit + 2)(Nu: + l)P(N H -. + 2,iV D: ) - N H: (N H: - l)P(N H: ,N D: )] 
+ A D [(iVb; + 2)(N D: + 1)P(N H ..,N D ; + 2) - N D ..(N D .. - l)P(N H , N D -.)] 
+ (A H + A D ) [{N h -. + l)(tf D: + l)P(N H -. + I, No-. + 1) - M ft P(iV K ,4)]. 



(36) 



The set of equations is written for the various integer values of A^h : and A^d- These equations are 
essentially identical to those derived previously by Lipshtat et al. ( 2004} , except for the fact that 
we have explicitely introduced the possibility of the Langmuir rejections via the —F ■ N/S terms. 
Using the relation /i = we rewrite the master equation in the following way: 



dP(N H ;,N D ;) 

dt 



F H [P(Nu: -l,N D: )- P(N H :,N D: )] 
+ F D [(P(Nh-,N D: - 1) - P(N H: ,N D ;] 
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+ fa [Nh P(N H :,N D ;) - (N H : - 1)P(N H ; - 1,N D ..)] 
+ / D [N H P(Nh:,N D ..) - N H .P(Nh:,N D .. - 1)] 

+ W H [(Nn: + 1)P(N H : + 1, AW - N H .P(N H .,,N D ..)] 

+ W D [(JVb: + l)P(Mi:,M: + 1) - N D .P(N H: ,N D: )] 

+ A H KM* + 2)(Nn + 1 )P(Nh: + 2,M>) - M*(M* - Wtt,W 

+ A D [(N D: + 2)(M: + l)P(N H ..,N D .. + 2) - M:(M>: - lWa,^)] 

+ (A H + A D ) [(M* + 1)(M, : + 1)P(M* + 1,M> + 1) - N H: N D: P(N H: ,N D: )] (37) 

This set of equations is identical to the standard master equations system ( |Barzel & Bih am 2007b) 
with two additional terms proportional to fn and fp. By a suitable summation over the probabil- 
ities in the master equation, one can obtain the moments of the distribution of adsorbed species 
populations, defined by: 

X W H:<:%:,iVD:) (38) 

Nh:,Nd; 

The order of each moment is defined by the sum / + k. The first-order moments (Mt : ) and (M>:) 
represent the mean number of adsorbed H and D atoms, respectively. 

3.2.2. Moment equations 



Lipshtat & Biham ( 2003 ) have shown that the moment equation formalism may describe adequately 



the evolution of the system and allows a significant reduction of the number of coupled equations 
that needs to be solved. We apply the same technique as in Barzel & Biham ( 2007b| l to derive the 



corresponding moment equations when the additional terms introduced by the Langmuir rejection 
are included. The time derivatives of the moments introduced by the term proportional to f is 



dt 

d(N D: ) f 
dt 



flKMt: 2 } 



2 M4:WH:,M:) = -/H<M4:> 
N h ,Nd 

N H ,N D 

Yj 4^H:,Wd:) = -M(Nk) + 2<M|» 

Nh,Nd 

2 N D .P f (N H ;,N D ..) = -/ D «Mi:> + 2{N H -N D :)) 

Nh,Nd 

Yj W„Mn:^/(Ml:,M>:) = _ /k(Ml:M):) - /d<< : >, (39) 
N H ,N D 

where the f dependent terms are given by 

Pf(N H :,N D ;) = / H [N H .P(Nh:,N D: ) - (Nh: ~ 1 Wh: - l.Afo:)] 

+ f D [N H P(Nh:,N D: ) - N H .P(Nh:,N D ; - 1)] (40) 



dt 

d{N D . 2 ) f 
dt 

d(N H: N D ..)f 
dt 



Adding these terms to the standard moment equations, one obtains the following system: 



12 
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d(N n } 
dt 

d(N D .} 

dt 
d(N H . 2 ) 
dt 

d{N D - 2 ) 
dt 

d(N H: N D ..) 
dt 



F H + (2A H - Wh - /h)<Wh:> - 2A ll (N 2 li .) - (A H + A D )(N H .N D ) 

F D - fuiNn) + (2A D - W D ){N D .) - 2A D << : > - (A H + A D )(N H: N D: ) 

F H + (2F H + W H + 4A h - /hX^Vh:) - (2W H + 4A H + 2/ H )<< : ) 

(A H + A D ){N H -N D :) 

F D - f D (N H :) + (2F D + W D + 4A D ){N D .) - (2W D + 4A D )<A£> 
(A H +A D + 2f D )(N H .N D .) 

A H + A D + /h)<A^H:A^D:) 



(41) 



This is the set of equations that is actually used in the interstellar chemistry model, in the case where 
rejection only due to H atoms is included. This system of coupled differential equations accounts 
for the populations of H and D atoms on a grain (given by the two first order moments) and for the 
reaction rates (given by the three second order moments). One then obtains the formation rates r^, 
thd and r£> 2 , of H2, HD and D2 , respectively, on a single grain (Lipshtat et al. 2004): 



r Hl = A H ((N H ) - (N H )) 
ruD = (Ah + A D )(N H N D ) 
r D2 = A d (<a£>-<JV d >). 



(42) 



The formation rate of H2, HD and D2 on grains is further obtained by summing over the number of 
grains. 

3.3. Master and moment equations with rejection terms for both H and D atoms 

When the Langmuir rejection terms due to both H and D atoms adsorbed on the surface are in- 
cluded, the master equation takes the form : 



dP(N H ,N D .) 
dt 



= F„ 



+ F r 



[l-f- ^)p(Nh,N d .. - 1) - (l - fi - f) P (NH,N D . } 
+ W H KNh: + 1)P(N H : + 1,N D ..) - N H -P(N H: ,N D: )] 

+ W D [(JVb: + 1)P(Nh:,N D .. + 1) - N D .P(Nh:,N D: )] 

+ A H M H: + 2)(N H ; + l)P(N H : + 2,N D: ) - N H: (N H .. - 1)P(N H , N D ..)] 
+ A D [(N D: + 2)(N D: + l)P(N n: ,N D ., + 2) - N D ..(N D: - l)P(iV H: ,iV D: )] 
+ (A H + A D ) [(Nh: + 1)(A^d : + VP(Nh: + hN D: + 1)-Nh N d .P(Nh:,N d ..)]. 



(43) 



Additional terms proportional to fn and fo are introduced in the master equations to account the 
Langmuir rejection due to D atoms. The derivation of the moment equations is very similar to the 
previously described procedure, leading to the following set of coupled equations for the first and 
second moments: 
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dt 
d(N D .} 

dt 
d(N H . 2 ) 
dt 

d{N D - 2 ) 
dt 

d(N H: N D ..) 
dt 



= F H + (2A H - Wh - fnWn.) - MN D: ) - 2A H <iV H .> - (A H + A D )(N H N D: ) 
= F D - f D {N H .) + (2 Ad - W D - f D ){N D .) - 2A D (N 2 .) - (A H + A D )(jV H: iV D: > 

= F H + (2F H + W H + 4A H - / H )<JV H :> - / H <Afo:> - (2W H + 4A H + 2/ H )<A&> 



- (A H +A D +2f H )(N H .N D .) 
= F D - f D (N H: ) + (2F D + W 

- (A H +A D +2/ D )<^ H :iVD:> 



- / D <Afo:> + (2^d + W D + 4A D - / d )<JVd:> " (2W D + 4A D + 2f D ){N D .) 



= F D < N H : > +F H {N D: ) - / D << : > - f H {N D .) - (W H + W D 
+ A H + A D + / H + /dWh-Nb:) 



This is the set of equations that is actually used in the interstellar chemistry model, in the 
case where rejection due to both H and D atoms is included. The production rates of H2, HD 
and D2 on a single grain are given by the same equations as those reported previously [Eq. (|42j] . 
Integrating over the grain size distribution, we obtain the formation rates (in cm _3 s _1 ) as well as 
the adsorption and desorption rates and the number of adsorbed hydrogen and deuterium atoms (in 
cm 4 ) as described in the next section. 

4. Calculations and results 

4. 1 . Fixed dust temperatures 

We consider first a diffuse cloud model whose parameters are given in Table [2] The chemical net- 
work comprises of 134 chemical species containing H, D, C, N, O, S and a typical metal "Fe" 
undergoing charge exchange reactions only in addition to photoionisation and electronic recombi- 
nation. 

We compare the results obtained via several approximations within the rate and the moment 
equations on the formation of H2, HD and D2 on amorphous carbon grains (see Table [T] for their 
properties) for different fixed relevant dust temperatures (between 8 and 20K) and their link to the 
gas phase chemistry. We use the Meudon PDR code ( Le Petit et al.|2 002, 2006) in which we have 



introduced the appropriate changes. The model corresponds to a slab of gas of 6.06 pc (correspond- 
ing to a total visual magnitude of 1) irradiated from both sides by the standard interstellar radiation 
field as given by Draine ( 1978] >. 



4.1 .1 . Incorporation of the moment equations into the Meudon PDR code 

We couple the full gas phase network and the surface chemistry of H2, HD and D2 in the PDR 



code by solving explicitely the steady state of Eq. 44 describing the evolution of the first and 
second moments of the distribution of adsorbed hydrogen and deuterium atoms on a single grain. 
The values of the first order moment give the average numbers of adsorbed H and D atoms on a 
single grain whereas no obvious physical meaning is associated to the second order moments. We 
perform the integration on the size distribution to get the appropriate quantities such as the number 
densities of adsorbed H, D and gas phase atomic and molecular number densities. First, we have 
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Table 2. The parameters used in the model. Numbers in parentheses refer to powers of ten. Njj is 
the total column density of protons (N(H) + 2 N(H2)) and Ay is the visual extinction. 



Gas parameter 


Dust parameter 


n n ( cm ^ } 


100 


n ■ (cw\^\ 
L *mm V" lll J 




T 

1 gas 


70 


n (cw\\ 








n 

H 


3.5 


X (Draine unit) 


1 


material 


amorphous carbon 


D/H 


1.5(-5) 


r 


1 


C/H 


1.34(-4) 


G 


0.01 


O/H 


3.19(-4) 


total visual magnitude 


1 


N/H 


7.5(-5) 


N H /A V 


1.871(21) cnr 2 /magnitude 


S/H 


1.86(-5) 


^dust 


fixed [8 - 20K] 


Fe/H 


1.5(-8) 







Table 3. The three models simulated using the rate equations and the moment equations and their 
legends. 



Model 


rate/moment 


No rejection 


Rejection for H only 


Rejection for H and D 


RA 


Rate equations 


yes 






RB 


Rate equations 




yes 




RC 


Rate equations 






yes 


MA 


Moment system 


yes 






MB 


Moment system 




yes 




MC 


Moment system 






yes 



checked the relevance of the treatment by considering only H2 formation by comparing with the 
analytical solution obtained from the resolution of the two coupled equations (Lipshtat & Biham 
[2003] l. 

4.1.2. Formation of H 2 and HD on grains 

We have examined the assumptions described above, using both the moment equations and the 
rate equations. More specifically, we compared three models: a model that includes no Langmuir 
rejection, a model that includes rejection only due to adsorbed H atoms and a model that includes 
rejection due to both H and D atoms adsorbed on grains. The three models are listed in table|3] 

We will present the results in the form of the standard rate law for the production of molecular 
hydrogen in interstellar clouds. In this rate law, the production rate Ru n (cm 3 s~') is expressed by 

R Hl = a(H 2 )«(H)« H , (44) 

where «h (cm 4 ) is the total density of H nuclei, in atomic and molecular forms. It can be ap- 
proximated by ;ih = «(H) + 2«(H2). The parameter aiHi) (cm 3 s~') is the effective rate co- 
efficient. In the case of diffuse clouds, this rate coefficient is often crudely approximated by 
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2 10 



■16 




dust 



Fig. 1. The formation rate a^T^) (cm 3 s _1 ) of H2 molecules on grains, obtained from the rate equa- 
tions (solid line) and from the moment equations (dashed lines), where the grain sizes follow the 
MRN size distribution. Results obtained from the three models of Table 3 are presented: Model A 
which includes no rejection (open circles), model B which includes rejection only due to adsorbed 
H atoms (open squares) and model C which includes rejection due to both adsorbed H and adsorbed 
D (open diamonds). 



a(H2) = 5x10 17 yT gas /300. The formation rate of HD molecules can be expressed in a simi- 
lar way, by 

Rhd = a(HD)n(D)« H , (45) 

The parameter a(HD) (cm 3 s -1 ) is the effective rate coefficient for HD formation on dust grains. 

In Fig. [T] we present the rate coefficient a(H2) for the formation of H2 molecules on grains and 
in Fig.|2]we present the rate coefficient a(HD) for the formation of HD molecules. These rate coef- 
ficients are plotted versus the assumed fixed temperature of the grains for the three models specified 
above, simulated by both the rate equations and the moment equations. The conditions used in the 
simulations are those obtained at the edge of the modelled cloud, i.e. Ay = 0. When the Langmuir 
rejection effect is taken into account, the formation of molecular hydrogen on grains is efficient 
within a narrow window of the grain temperatures, in agreement with previous studies (Katz et al. 
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Fig. 2. The formation rate ff(HD) (cm 3 s _1 ) of HD molecules on grains, obtained from the rate 
equations and the moment equations, where the grain sizes follow the MRN size distribution. Same 
conventions as in Fig.[T] 



1999 Lipshtat et al.|2004) >. At higher temperatures, the hydrogen atoms do not stay long enough 



on the grain to encounter each other and form molecules. At lower temperatures, the grain surface 
is saturated by nearly immobile hydrogen atoms, and the reaction rate is very low. These adsorbed 
atoms reject new atoms that impinge on the surface and prevent the formation of additional layers 
of hydrogen atoms on the surface. When the Langmuir rejection term is removed from the equa- 
tions, the range of high efficient recombination is extended towards lower temperatures. However, 
in this case the model enables the accumulation of many layers of adsorbed hydrogen atoms on the 
surface. Such accumulation is unphysical. Furthermore, under such circumstances, the formation of 
molecular hydrogen is not mediated by diffusion, and thus the model itself is not valid. Laboratory 
experiments provide strong evidence for Langmuir rejection ( |Katz et al.|19 99 ). Therefore, the stan- 
dard forms of the rate equations and of the moment equations which we incorporate in the PDR 
code are those that include the Langmuir rejection term. In presence of this term, we observe that 
the high efficiency window is somewhat larger for HD formation than for H2 formation. This is 
due to the somewhat higher desorption barrier for D atoms than for H atoms. The maximum value 
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Table 4. Column densities (in cm 2 ) of H, H2, D and HD and the molecular fraction f for model 
C, obtained from the rate equations (RC) and from the moment equations (MC) for different grain 
temperatures. Numbers in parentheses refer to powers of ten. 



Model 


Tdust 


H 


Hi 


f 


D 


HD 


RC 


8 


1.87(21) 


1.60(14) 


1.7(-7) 


2.76(16) 


2.64(10) 


MC 


8 


1.87(21) 


1.60(14) 


1.7(-7) 


2.81(16) 


6.26(9) 


RC 


10 


1.77(21) 


5.30(19) 


5.65(-2) 


1.74(16) 


1.07(16) 


MC 


10 


1.84(21) 


1.62(19) 


1.73(-2) 


2.28(16) 


5.25(15) 


RC 


12 


6.65(19) 


9.02(20) 


0.96 


1.26(16) 


1.55(16) 


MC 


12 


6.33(19) 


9.04(20) 


0.97 


1.26(16) 


1.55(16) 


RC 


14 


6.29(19) 


9.04e(20) 


0.97 


1.26(16) 


1.55(16) 


MC 


14 


6.29(19) 


9.04(20) 


0.97 


1.26(16) 


1.55(16) 


RC 


16 


6.95(19) 


9.01(20) 


0.96 


1.25(16) 


1.56(16) 


MC 


16 


8.87(19) 


8.91(20) 


0.95 


1.21(16) 


1.60(16) 


RC 


18 


1.76(20) 


8.48(20) 


0.91 


1.07(16) 


1.74(16) 


MC 


18 


3.43(20) 


7.64(20) 


0.82 


8.91(15 


1.91(16) 


RC 


20 


9.60(20) 


4.55(20) 


0.49 


7.58(15) 


2.05(16) 


MC 


20 


1.52(21) 


1.77(20) 


0.19 


1.12(16) 


1.69(16) 



of the formation rate of H2shown in Fig. 1 is in excellent agreement with the "effective" formation 
rate osh 2 (cm 3 s _1 ), expressed by (Le Bourlot et al. 



1995) 



a Hl = 1.4 3 • m H G ^ = ^ ^ ^-16^ (46) 

4 • p yj^min^max 

with the chosen parameters. The corresponding value for HD is ohd - VSothj = 2.3 x 10 _16 «h 
(cm 3 s _1 ), as seen on Fig. [2] The rate equations and the moment equations are almost identical at the 
maximum when desorption effects are not significant. However, significant discrepancies appear at 
the edges. Rejection effects play a significant role at the lowest temperatures of the grains, which 
may not be very relevant for astrophysical purposes. We display in Table |4] the resulting column 
densities of gas phase H, H2, D and HD when integration of the abundance densities is performed 
over the width of the cloud. We have chosen to display only the results of model C which are the 
most accurate as rejection effects are included both for H and D impinging atoms. 

The column densities of H2 display significant variations with the assumed dust temperature as 
reflected by the values of the molecular fraction defined by / = 2N(H 2 )/[N(H) + 2N(H 2 )]. The 
resulting column densities of HD become less sensitive as soon as some molecular hydrogen is 
present. This is due to the fact that in presence of H2 the formation of HD occurs mainly via gas 
phase processes. The formation of HD on grains is more efficient than in the gas phase only at 
the edge where no H2 is yet formed. We discuss this in more detail in the next section where we 
introduce a variable dust temperature as a function of the visual magnitude. 
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4.2. A v dependent dust temperatures 



4.2.1 . Homogeneous temperature distribution 

It is well recognised that the grain temperatures decrease when the strength of the impinging ra- 
diation field decreases for increasing visual magnitude. We introduce such a dependence by using 



the simple analytic formula presented by Hollenbach et al. ( 1991 1, which provides the grain tern 



perature as a function of the strength of the radiation field and the visual magnitude. As our model 
involves a slab of gas, irradiated from both sides, we extend the formula of Hollenbach et al. ( 1991) ) 
by introducing \xieft ex P( _ 1 -8A V ) +Xog& ex P{- 1 .8(Aj, ot - A v )] for the A v dependence, where xieft and 
fright are the strengths of the radiation field on both sides, expressed in Draine units, and A™ is the 
total visual magnitude of the cloud. In this approximation, grains of all sizes exhibit the same 



average temperature. We consider possible fluctuations effects in 4.2.2 
Table 5. Physical conditions relevant to a dense PDR. 



Physical quantity value 

n H 10 4 crrT 3 

Af 10 

size 0.61 pc 

ATleft 100 

/bright 1 

T m 30K 

Tdust Ay dependent 
size distribution of grains same as in Tablej^J 



In order to span a significant range of dust temperatures, we consider a dense photodissociation 
region (PDR) with a radiation field strength (xieft) of 100 impinging on the left side and the standard 
radiation field Cfright) impinging on the right side. The parameters of the model are displayed in 
Table [5] As we want here to discuss the different effects of the rate versus moment equations 
systems, we keep all the physical parameters identical. 

We display in Figures[3]and|4]the dust temperature variation together with the formation rates 
through gas phase and grain surface reactions of H2 and HD obtained with the RC and MC treat- 
ments as these are expected to be the most physical within the rate equation and moment equation 
formalisms. The dust temperature profile spans a range of values between about 27K on the left 
side, reaching a minimum of 10K in the shielded region at Ay between 5-7 , and reaching a value of 
about 10K on the right edge. This range corresponds to very low to high efficiencies of formation 
of H2 and HD via grain surface reactions; As far as H2 formation is concerned, we see that the 
difference between rate and moment equation formalisms is significant at the left side of the cloud 
where the dust temperature reaches a value of 27K, corresponding to a situation where H atoms 
desorb efficiently from the surface. The results obtained via the rate equations are larger by a factor 
of about 2. In this particular case, this leads to a competition between gas phase and dust processes. 
However, as soon as the dust temperature reaches a value of about 25, the formation efficiency on 
dust becomes preponderant. Both treatments are equivalent for a range of temperature between 20 
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Fig. 3. Comparison between gas phase formation (circles) and grain-surface formation (squares) 
rates of H2 (in cm _3 s _1 ), obtained from the rate equations (solid lines and full symbols) and from 
the moment equations (dotted lines and open symbols). The thin solid line shows the grain temper- 
ature (the temperature scale is on the right hand side) 



Table 6. Column densities in cm 2 and molecular fraction f for the dense PDR models. Numbers 
in parentheses refer to power of ten. 



Model 


H 


H; 


f 


D 


HD 


C 


CO 


CH 


RC 


2.29(21) 


8.21(21) 


0.88 


3.29(16) 


2.28(17) 


7.74(17) 


1.18(18) 


9.31(14) 


MC 


6.35(21) 


6.18(21) 


0.66 


3.25(16) 


2.46(17) 


5.1(17) 


1.10(18) 


1.02(15) 



and 1 IK as already shown on Figure[T| Significant discrepancies occur again for the range of visual 
magnitude 5-8 where the dust temperature reaches a value of 10K. 

These behaviours have a direct impact on the HD formation mechanisms, as displayed on 
Figure [4] As long as no significant formation of H2 occurs, HD is mainly formed via dust pro- 
cesses. However, as soon as H2 is present, formation of HD via gas phase is much more efficient. 
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Fig. 4. Comparison between gas phase formation (circles) and grain-surface formation (squares) 
rates of HD (in cm ^s 1 ) obtained from the rate equations (solid lines and full symbols) and from 
the moment equations (dotted lines and open symbols). The thin solid line shows the grain temper- 
ature (the temperature scale is on the right hand side) 



In the shielded region where Ay reaches a value of a few and the dust temperature is about 10K, 
significant variations occur due to rejection effects. Table|6]displays the resulting column densities. 
Whereas significant differences arise for the column densities of H and H2, the values relevant to 
the other species are quite similar. 

4.2.2. Temperature fluctuations effects 

The consideration of a homogeneous distribution of dust temperatures as a function of Ay is rea- 



sonable for low impinging radiation fields and grains of radii larger than 1 .0 x 10 6 cm, ( Horn et al. 
2007 ). For grains of radii smaller than ~ 1.0 x 10~ 6 cm, the absorption of a single UV photon may 



cause a temperature spike that heats up the grains by over 20 K ( Horn et al.|2007 1. Such a spike is 



likely to lead to the immediate desorption of all the H and D atoms adsorbed on the grain surface. 
This conclusion is also derived from detailed kinetic Monte Carlo simulations of H2 formation on 
stochastically heated olivine grains (Cuppe rTet al.|200 6; Herbst & Cuppen 2006). To quantify the 
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effect of temperature fluctuations on small grains, we have run similar models than those presented 
previously, except that we cut the size distribution of grains at the minimum value of 1.0 x 10~ 6 
cm. We present in Table [7] the corresponding column densities. In this way, we do not take into 
account any production of molecular hydrogen on grains of radii smaller than 1.0 x 10 6 cm both 
in rate equations and systems of moments approach. In these models, computed column densities 

Table 7. Column densities in cirT 2 and molecular fraction f for the dense PDR models in which 
surface chemistry is cut for very small grains (r < 1.0 x 10~ 6 cm). Numbers in parentheses refer to 
power of ten. 

Model H H 2 f D HD C CO CH 

RC 5.88(21) 6.41(21) 0.69 3.52(16) 2.43(17) 8.45(17) 1.11(18) 1.02(15) 
MC 8.40(21) 5.16(21) 0.55 3.47(16) 2.46(17) 9.27(17) 1.02(18) 1.06(15) 

of H2 and HD are smaller than those obtained with the MRN distribution extended to 3.0 x 10~ 7 
cm (cf Table[7]i due to a smaller integrated grain cross section value. We find that rate and moment 
equations give still significantly different values of the total column densities of H2 even for these 
larger grains. These differences arise essentially from the regions when T^, is larger than about 
18K, i.e. at the edges of the cloud. 

The present results give an upper limit of the effect of temperature fluctuations as these can only 
result in a decrease of molecule formation at the surface of grains. In addition, these effects will not 
occur at A v values larger than ~ 1 and smaller than ~ 9 as the radiation field is then considerably 
reduced. 

5. Summary 

We have incorporated the moment equations for the formation of H2 and its deuterated versions 
into the Meudon PDR code, examined their applicability and relevance and compared the results 
with those obtained from rate equations incorporated into the same code under identical conditions. 
Previous analysis has shown that as long as the temperatures of the dust grains are not too high, the 
reaction rates obtained from the moment equations coincide with those obtained from the rate equa- 
tions. At grain temperatures around 18 K or higher, on the high-temperature end of the efficiency 
window there are significant differences, where the rate equations over-estimate the formation rates 
of H2 and HD. These deviations are mainly due to the very small grains, on which the population 
sizes of adsorbed H and D atoms exhibit large fluctuations. Under such conditions it is important to 
use the moment equations rather than the rate equations. For low grain temperatures, below 12 K or 
so, we find that the Langmuir rejection term makes a crucial difference. In this range, the rejection 
term is required in order to prevent the freezing of multi-layers of hydrogen atoms on the grains. 
Such freezing is unphysical. Here we have derived for the first time a set of moment equations in 
which the Langmuir rejection term is incorporated. The comparison was made with rate equations, 
in which a suitable rejection term was also incorporated. For very low temperatures, where the 
grains are covered by a layer of H and D atoms and the rejection term plays a major role, we also 
find a discrepancy between the formation rates of H2 and HD obtained from the rate equations 
and the moment equations. Under these conditions, deep inside a molecular cloud, the formation 
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rate and dissociation rate of H2 molecules are low, while HD molecules form mainly by gas phase 
reactions. These are conditions under which both the rate equations and the moment equations are 
accurate. The observed discrepancy is due to somewhat different physical consitions that emerge 
as a result of the discrepancy in warmer regions, where the moment equations are accurate and the 
rate equations are not. Therefore, the moment equation results are those one can rely on also in the 
cold regions. 

We conclude that the moment equations provide an efficient method for the evaluation of molec- 
ular hydrogen formation in interstellar chemistry codes. The equations provide accurate results for 
the reaction rates on both large and small grains. The equations are easy to construct and are ef- 
ficient in terms of computational resources. They can be easily coupled to the rate equations of 
gas-phase chemistry. Since the moment equations are linear, their stability and convergence prop- 
erties are often superior even in comparison to the rate equations. 

The present study has been achieved by considering that molecular hydrogen can be formed 
only through diffusion of adsorbed atoms, with the strong constraints of present experimental 



knowledge. Other formation mechanisms may be at work, as those suggested by Cazaux & Tielens 



( 2004 ) and Habart et al. ( 2004 1 in their interpretation of observations of warm molecular hydrogen. 



The method can now be extended to more complex reaction networks on grains. The network 
that involves H, O and CO molecules, from which ice mantles containing H2O, CO2 and CH3OH 



are formed by successive hydrogenation and oxydation reactions, is promising (Barzel & Biham 



2007b I. Unlike the H and D atoms, these heavier atoms and molecules are more strongly bound 
to the surface and to each other and do not exhibit the Langmuir rejection property. Therefore, in 
the construction of moment equations for more complex networks one only needs to maintain the 
rejection terms for H and D, presented above and there is no need to incorporate such terms for any 
other atomic or molecular species. Recent studies have shown that CH3OH molecules do not form 
in the gas phase as efficiently as previously expected. The approach presented here will enable to 
evaluate the formation rate of CH3OH on dust grains under different physical conditions, taking 
into account the contributions of grains of different sizes. 
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